Mapping guides: -More general intro: https://www.joshuastevens.net/cartography/make-a-bivariate-choropleth-map/ -R tutorial: https://cran.r-project.org/web/packages/biscale/vignettes/biscale.html
library(tidyverse)
library(leaflet)
library(sf)
## Linking to GEOS 3.8.1, GDAL 3.1.1, PROJ 6.3.1
library(rgdal)
## Loading required package: sp
## rgdal: version: 1.5-18, (SVN revision 1082)
## Geospatial Data Abstraction Library extensions to R successfully loaded
## Loaded GDAL runtime: GDAL 3.1.1, released 2020/06/22
## Path to GDAL shared files: /Library/Frameworks/R.framework/Versions/4.0/Resources/library/rgdal/gdal
## GDAL binary built with GEOS: TRUE
## Loaded PROJ runtime: Rel. 6.3.1, February 10th, 2020, [PJ_VERSION: 631]
## Path to PROJ shared files: /Library/Frameworks/R.framework/Versions/4.0/Resources/library/rgdal/proj
## Linking to sp version:1.4-4
## To mute warnings of possible GDAL/OSR exportToProj4() degradation,
## use options("rgdal_show_exportToProj4_warnings"="none") before loading rgdal.
# read in project data
merge_data = read_csv("data/merge_data.csv") %>%
rename(STUSPS = state_id) %>% # to match shapefile col name below
mutate_if(is.numeric, round, digits = 2) # round to 2 decimal places
##
## ── Column specification ────────────────────────────────────────────────────────
## cols(
## .default = col_double(),
## state_id = col_character(),
## hostility = col_character()
## )
## ℹ Use `spec()` for the full column specifications.
merge_data
## # A tibble: 50 x 27
## STUSPS counties_no_pro… percent_abortion percent_birth percent_medicaid
## <chr> <dbl> <dbl> <dbl> <dbl>
## 1 AK 79 11 73 24
## 2 AL 91 10 74 18
## 3 AR 97 9 75 28
## 4 AZ 80 11 73 24
## 5 CA 24 20 65 27
## 6 CO 77 13 71 20
## 7 CT 13 23 62 23
## 8 DE 33 19 66 19
## 9 FL 67 21 64 17
## 10 GA 94 17 68 13
## # … with 40 more rows, and 22 more variables: percent_women_no_provider <dbl>,
## # percent_uninsured <dbl>, percent_bc_18_49 <dbl>,
## # abortion_rate_15_17_state <dbl>, abortion_rate_15_19_state <dbl>,
## # abortion_rate_15_44_state <dbl>, abortion_rate_18_19_state <dbl>,
## # birthrate_15_17_state <dbl>, birthrate_15_19_state <dbl>,
## # birthrate_18_19_state <dbl>, need_bc_hisp_20_44 <dbl>,
## # need_bc_hisp_younger_20 <dbl>, need_bc_13_44 <dbl>,
## # need_bc_black_20_44 <dbl>, need_bc_black_under_20 <dbl>,
## # need_bc_white_20_44 <dbl>, need_bc_white_under_20 <dbl>,
## # need_bc_under_20 <dbl>, total_expend_abortion <dbl>,
## # public_expenditures <dbl>, expenditure_rate <dbl>, hostility <chr>
# import state shapefile from data folder
all_states <- readOGR("data/cb_2018_us_state_500k/cb_2018_us_state_500k.shp")
## OGR data source with driver: ESRI Shapefile
## Source: "/Users/rht2112/Documents/My Documents/Grad School/Coursework/Semester 3/Data Science I/Homeworks and Projects/p8105_2020_final/data/cb_2018_us_state_500k/cb_2018_us_state_500k.shp", layer: "cb_2018_us_state_500k"
## with 56 features
## It has 9 fields
## Integer64 fields read as strings: ALAND AWATER
summary(all_states)
## Object of class SpatialPolygonsDataFrame
## Coordinates:
## min max
## x -179.1489 179.77847
## y -14.5487 71.36516
## Is projected: FALSE
## proj4string : [+proj=longlat +datum=NAD83 +no_defs]
## Data attributes:
## STATEFP STATENS AFFGEOID GEOID
## Length:56 Length:56 Length:56 Length:56
## Class :character Class :character Class :character Class :character
## Mode :character Mode :character Mode :character Mode :character
## STUSPS NAME LSAD ALAND
## Length:56 Length:56 Length:56 Length:56
## Class :character Class :character Class :character Class :character
## Mode :character Mode :character Mode :character Mode :character
## AWATER
## Length:56
## Class :character
## Mode :character
data_sf <- merge(all_states, merge_data, all.x = F)
summary(data_sf)
## Object of class SpatialPolygonsDataFrame
## Coordinates:
## min max
## x -179.14891 179.77847
## y 18.91036 71.36516
## Is projected: FALSE
## proj4string : [+proj=longlat +datum=NAD83 +no_defs]
## Data attributes:
## STUSPS STATEFP STATENS AFFGEOID
## Length:50 Length:50 Length:50 Length:50
## Class :character Class :character Class :character Class :character
## Mode :character Mode :character Mode :character Mode :character
##
##
##
##
## GEOID NAME LSAD ALAND
## Length:50 Length:50 Length:50 Length:50
## Class :character Class :character Class :character Class :character
## Mode :character Mode :character Mode :character Mode :character
##
##
##
##
## AWATER counties_no_provider percent_abortion percent_birth
## Length:50 Min. :13.00 Min. : 5.00 Min. :57.00
## Class :character 1st Qu.:75.50 1st Qu.: 9.00 1st Qu.:68.00
## Mode :character Median :91.00 Median :11.50 Median :72.50
## Mean :78.24 Mean :13.36 Mean :71.12
## 3rd Qu.:95.00 3rd Qu.:17.00 3rd Qu.:75.00
## Max. :99.00 Max. :29.00 Max. :79.00
##
## percent_medicaid percent_women_no_provider percent_uninsured percent_bc_18_49
## Min. : 9.00 Min. : 0.00 Min. : 3.00 Min. :61.80
## 1st Qu.:15.50 1st Qu.:18.25 1st Qu.: 7.00 1st Qu.:67.90
## Median :19.50 Median :42.50 Median :11.00 Median :71.10
## Mean :20.24 Mean :43.16 Mean :10.82 Mean :70.44
## 3rd Qu.:24.75 3rd Qu.:62.75 3rd Qu.:13.00 3rd Qu.:73.45
## Max. :38.00 Max. :96.00 Max. :24.00 Max. :78.00
## NA's :12
## abortion_rate_15_17_state abortion_rate_15_19_state abortion_rate_15_44_state
## Min. :1.000 Min. : 2.600 Min. : 4.600
## 1st Qu.:2.200 1st Qu.: 4.600 1st Qu.: 7.675
## Median :3.050 Median : 6.000 Median : 9.900
## Mean :3.468 Mean : 6.684 Mean :11.238
## 3rd Qu.:4.450 3rd Qu.: 8.200 3rd Qu.:14.100
## Max. :9.600 Max. :17.200 Max. :28.200
##
## abortion_rate_18_19_state birthrate_15_17_state birthrate_15_19_state
## Min. : 4.80 Min. : 3.500 Min. : 8.50
## 1st Qu.: 8.20 1st Qu.: 6.400 1st Qu.:15.93
## Median :10.00 Median : 8.250 Median :19.75
## Mean :11.49 Mean : 8.648 Mean :20.72
## 3rd Qu.:13.25 3rd Qu.:10.175 3rd Qu.:24.80
## Max. :28.40 Max. :15.000 Max. :34.70
##
## birthrate_18_19_state need_bc_hisp_20_44 need_bc_hisp_younger_20
## Min. :13.80 Min. : 680 Min. : 270
## 1st Qu.:29.75 1st Qu.: 11792 1st Qu.: 2745
## Median :37.85 Median : 28845 Median : 7275
## Mean :38.99 Mean : 82409 Mean : 19617
## 3rd Qu.:48.75 3rd Qu.: 66735 3rd Qu.: 16162
## Max. :66.50 Max. :1102300 Max. :263760
##
## need_bc_13_44 need_bc_black_20_44 need_bc_black_under_20
## Min. : 34960 Min. : 160 Min. : 60
## 1st Qu.: 119260 1st Qu.: 4402 1st Qu.: 1148
## Median : 303910 Median : 26175 Median : 6650
## Mean : 411957 Mean : 59597 Mean :14382
## 3rd Qu.: 463600 3rd Qu.:101400 3rd Qu.:24575
## Max. :2526010 Max. :264700 Max. :61270
##
## need_bc_white_20_44 need_bc_white_under_20 need_bc_under_20
## Min. : 11630 Min. : 2310 Min. : 8220
## 1st Qu.: 52142 1st Qu.: 18342 1st Qu.: 26780
## Median :121200 Median : 41105 Median : 68105
## Mean :148834 Mean : 52164 Mean : 92522
## 3rd Qu.:213080 3rd Qu.: 71735 3rd Qu.:103850
## Max. :468270 Max. :161870 Max. :527440
##
## total_expend_abortion public_expenditures expenditure_rate hostility
## Min. : 0 Min. : 762 Min. :0.0100 Length:50
## 1st Qu.: 0 1st Qu.: 8054 1st Qu.:0.0500 Class :character
## Median : 76 Median : 24971 Median :0.0900 Mode :character
## Mean : 1931 Mean : 41698 Mean :0.0924
## 3rd Qu.: 453 3rd Qu.: 56273 3rd Qu.:0.1200
## Max. :32613 Max. :454706 Max. :0.2700
## NA's :13
# make color palette
abortion_rate_pal <- colorNumeric(
palette = "inferno",
domain = data_sf$percent_abortion,
reverse = TRUE)
# base map
base_map <- leaflet(data_sf) %>%
addProviderTiles("CartoDB.PositronNoLabels")
# add in custom labels for percent abortion
labels_abortion_rate <- sprintf(
"<strong>%s</strong><br/>%g%",
data_sf$STUSPS, data_sf$percent_abortion
) %>% lapply(htmltools::HTML)
# full map
abortion_map = base_map %>%
addPolygons(
fillColor = ~abortion_rate_pal(data_sf$percent_abortion),
weight = 2,
opacity = 1,
color = "white",
dashArray = "3",
fillOpacity = 0.7,
highlight = highlightOptions(
weight = 5,
color = "#666",
dashArray = "",
fillOpacity = 0.7,
bringToFront = TRUE),
label = labels_abortion_rate,
labelOptions = labelOptions(
style = list("font-weight" = "normal", padding = "3px 8px"),
textsize = "15px",
direction = "auto")) %>%
addLegend(pal = abortion_rate_pal, values = ~percent_abortion, opacity = 0.7, title = 'Abortion Rates by State <br> (% of pregnancies ending in abortions)',
position = "bottomright")
abortion_map
# make color palette
bc_fund_pal <- colorNumeric(
palette = "inferno",
domain = data_sf$expenditure_rate,
reverse = TRUE)
# add in custom labels for public contraceptive expenditure
labels_bc_fund <- sprintf(
"<strong>%s</strong><br/>$%g/woman",
data_sf$STUSPS, data_sf$expenditure_rate
) %>% lapply(htmltools::HTML)
# full map
bc_fund_map = base_map %>%
addPolygons(
fillColor = ~bc_fund_pal(data_sf$expenditure_rate),
weight = 2,
opacity = 1,
color = "white",
dashArray = "3",
fillOpacity = 0.7,
highlight = highlightOptions(
weight = 5,
color = "#666",
dashArray = "",
fillOpacity = 0.7,
bringToFront = TRUE),
label = labels_bc_fund,
labelOptions = labelOptions(
style = list("font-weight" = "normal", padding = "3px 8px"),
textsize = "15px",
direction = "auto")) %>%
addLegend(pal = bc_fund_pal, values = ~expenditure_rate, opacity = 0.7, title = 'Public Contraceptive Expenditure Rate by State <br> ($/woman in likely need of publicly-funded services)',
position = "bottomright")
bc_fund_map
# make color palette
abortion_access_pal <- colorNumeric(
palette = "inferno",
domain = data_sf$percent_women_no_provider,
reverse = TRUE)
# add in custom labels for public contraceptive expenditure
labels_abortion_access <- sprintf(
"<strong>%s</strong><br/>%g%",
data_sf$STUSPS, data_sf$percent_women_no_provider
) %>% lapply(htmltools::HTML)
# full map
abortion_access_map = base_map %>%
addPolygons(
fillColor = ~abortion_access_pal(data_sf$percent_women_no_provider),
weight = 2,
opacity = 1,
color = "white",
dashArray = "3",
fillOpacity = 0.7,
highlight = highlightOptions(
weight = 5,
color = "#666",
dashArray = "",
fillOpacity = 0.7,
bringToFront = TRUE),
label = labels_abortion_access,
labelOptions = labelOptions(
style = list("font-weight" = "normal", padding = "3px 8px"),
textsize = "15px",
direction = "auto")) %>%
addLegend(pal = abortion_access_pal, values = ~percent_women_no_provider, opacity = 0.7, title = 'Percent of women living in a county <br> without access to abortion clinic by State',
position = "bottomright")
abortion_access_map
# make color palette
teen_births_pal <- colorNumeric(
palette = "inferno",
domain = data_sf$birthrate_15_19_state,
reverse = TRUE)
# add in custom labels for public contraceptive expenditure
labels_teen_births <- sprintf(
"<strong>%s</strong><br/> %g births per 1000 women",
data_sf$STUSPS, data_sf$birthrate_15_19_state
) %>% lapply(htmltools::HTML)
# full map
teen_births_map = base_map %>%
addPolygons(
fillColor = ~teen_births_pal(data_sf$birthrate_15_19_state),
weight = 2,
opacity = 1,
color = "white",
dashArray = "3",
fillOpacity = 0.7,
highlight = highlightOptions(
weight = 5,
color = "#666",
dashArray = "",
fillOpacity = 0.7,
bringToFront = TRUE),
label = labels_teen_births,
labelOptions = labelOptions(
style = list("font-weight" = "normal", padding = "3px 8px"),
textsize = "15px",
direction = "auto")) %>%
addLegend(pal = teen_births_pal, values = ~birthrate_15_19_state, opacity = 0.7, title = 'No. of births per 1000 women aged 15-19',
position = "bottomright")
teen_births_map
# make color palette
need_bc_teen_pal <- colorNumeric(
palette = "inferno",
domain = data_sf$need_bc_under_20,
reverse = TRUE)
# add in custom labels for public contraceptive expenditure
labels_need_bc_teen <- sprintf(
"<strong>%s</strong><br/>%g",
data_sf$STUSPS, data_sf$need_bc_under_20
) %>% lapply(htmltools::HTML)
# full map
need_bc_teen_map = base_map %>%
addPolygons(
fillColor = ~need_bc_teen_pal(data_sf$need_bc_under_20),
weight = 2,
opacity = 1,
color = "white",
dashArray = "3",
fillOpacity = 0.7,
highlight = highlightOptions(
weight = 5,
color = "#666",
dashArray = "",
fillOpacity = 0.7,
bringToFront = TRUE),
label = labels_need_bc_teen,
labelOptions = labelOptions(
style = list("font-weight" = "normal", padding = "3px 8px"),
textsize = "15px",
direction = "auto")) %>%
addLegend(pal = need_bc_teen_pal, values = ~need_bc_under_20, opacity = 0.7, title = 'No. of women who likely need public support for contraceptive services and supplies younger than 20',
position = "bottomright")
need_bc_teen_map
# make color palette
need_bc_pal <- colorNumeric(
palette = "inferno",
domain = data_sf$need_bc_13_44,
reverse = TRUE)
# add in custom labels for public contraceptive expenditure
labels_need_bc <- sprintf(
"<strong>%s</strong><br/>%g",
data_sf$STUSPS, data_sf$need_bc_13_44
) %>% lapply(htmltools::HTML)
# full map
need_bc_map = base_map %>%
addPolygons(
fillColor = ~need_bc_pal(data_sf$need_bc_13_44),
weight = 2,
opacity = 1,
color = "white",
dashArray = "3",
fillOpacity = 0.7,
highlight = highlightOptions(
weight = 5,
color = "#666",
dashArray = "",
fillOpacity = 0.7,
bringToFront = TRUE),
label = labels_need_bc,
labelOptions = labelOptions(
style = list("font-weight" = "normal", padding = "3px 8px"),
textsize = "15px",
direction = "auto")) %>%
addLegend(pal = need_bc_pal, values = ~need_bc_13_44, opacity = 0.7, title = 'Total No. of women who likely need public support for contraceptive services and supplies',
position = "bottomright")
need_bc_map
# make color palette
expenditure_rate_pal <- colorNumeric(
palette = "inferno",
domain = data_sf$expenditure_rate,
reverse = TRUE)
# add in custom labels for public contraceptive expenditure
labels_expenditure_rate <- sprintf(
"<strong>%s</strong><br/>%g",
data_sf$STUSPS, data_sf$expenditure_rate
) %>% lapply(htmltools::HTML)
# full map
expenditure_rate_map = base_map %>%
addPolygons(
fillColor = ~expenditure_rate_pal(data_sf$expenditure_rate),
weight = 2,
opacity = 1,
color = "white",
dashArray = "3",
fillOpacity = 0.7,
highlight = highlightOptions(
weight = 5,
color = "#666",
dashArray = "",
fillOpacity = 0.7,
bringToFront = TRUE),
label = labels_expenditure_rate,
labelOptions = labelOptions(
style = list("font-weight" = "normal", padding = "3px 8px"),
textsize = "15px",
direction = "auto")) %>%
addLegend(pal = expenditure_rate_pal, values = ~expenditure_rate, opacity = 0.7, title = 'Ratio of Total reported public expenditures for family planning client services in 1000s of dollars to number of women in likely need of public support for contraceptive services and supplies',
position = "bottomright")
expenditure_rate_map